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Abstract 



Pattern formation in parametric surface waves is studied in the limit of 
weak viscous dissipation. A set of quasi-potential equations (QPEs) is intro- 
duced that admits a closed representation in terms of surface variables alone. 
A multiscale expansion of the QPEs reveals the importance of triad reso- 
nant interactions, and the saturating effect of the driving force leading to a 
gradient amplitude equation. Minimization of the associated Lyapunov func- 
tion yields standing wave patterns of square symmetry for capillary waves, 
and hexagonal patterns and a sequence of quasi-patterns for mixed capillary- 
gravity waves. Numerical integration of the QPEs reveals a quasi-pattern of 
eight-fold symmetry in the range of parameters predicted by the multiscale 
expansion. 
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When a fluid layer with a free upper surface is subjected to vertical oscillation, Faraday 
waves are observed ]I]|§. In a large enough container and for low viscosity fluids, standing 
wave patterns of square symmetry are observed near threshold ||. Based on amplitude 
equations that we derive below, the experimental observation of square patterns in the 
capillary-dominated regime is explained. We also predict that hexagonal and a sequence 
of quasi-patterns can be stabilized for the case of a sinusoidal driving force as a result of 
triad resonant interactions for mixed capillary-gravity waves @. To our knowledge, this is 
the first theoretical derivation that starting from a realistic model of the fluid shows that a 
quasi-crystalline pattern is a stable steady state, and corroborates the conjecture of Newell 
and Pomeau || on the existence of the so-called "turbulent crystals" . 

Pattern-forming instabilities occur in a variety of extended nonlinear systems. The emer- 
gence of spatial patterns close to onset of the instability can often be described by amplitude 
equations However, for near-Hamiltonian (or weakly dissipative) systems, there is no 

general agreement on how dissipation should be incorporated into the amplitude equation 
formulation. Previous work on Faraday waves was based on amplitude equations for a purely 
Hamiltonian system, to which linear and nonlinear damping terms were added by introducing 
a dissipation function |7||| . In this approach, linear dissipative effects in the original system 
contribute only to linear damping terms in the associated amplitude equations, while nonlin- 
ear damping terms result entirely from nonlinear dissipative effects. Such an approach has 
contributed to the general belief that for near-Hamiltonian systems, nonlinear saturation of 
the linear instability does not occur if only linear dissipative effects are considered, and weak 
nonlinear dissipative or other higher order effects are needed for nonlinear saturation 
In this paper, we show that in the case of weakly damped parametric surface waves linear 
dissipative effects do contribute to the nonlinear damping terms in the amplitude equation, 
and that they alone can saturate the parametric instability. In addition, the experimental 
observation of square patterns in capillary-dominated regime is naturally explained without 
having to invoke poorly understood nonlinear dissipative effects, or higher order terms in 
the amplitude equation. 
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The basic difference with previous studies [0,|J is that although the bulk flow does re- 
main potential, it is modified by a rotational viscous boundary layer near the free surface 
that has to be explicitly incorporated into the analysis fT0j-|T2||. When the thickness of the 
viscous boundary layer is small compared to the typical wavelength of the pattern, the weak 
effects due to viscosity can be taken into account by introducing effective boundary con- 
ditions for the otherwise potential bulk flow. This is the basic idea of the quasi-potential 
approximation introduced below. We first expand the equations governing the motion of 
an incompressible viscous fluid and the appropriate boundary conditions at the free surface 
in the small thickness of the free surface boundary layer, 5. The resulting equations are 
further simplified by recasting them in a nonlocal form that involves the flow variables on 
the free surface only; thus eliminating the need to explicitly solve for the flow in the bulk 



1"2|| . Let z be the normal direction to the surface at rest and g(t) = —g — g z0 smttt the 
driving force where go is the constant acceleration of gravity, and Q and g z o are the angular 
frequency and the amplitude of the driving force respectively. We choose 1/uj = 2/fl as the 
unit of time and l/k as the unit of length, with k Q defined by uJq = gok + Tk^/ p, the linear 
dispersion relation for surface waves, where T is the surface tension and p the density of the 
fluid. Then the dimensionless, nonlocal and quasi-potential equations read JT2[ , 
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-(£>$) 2 -(£>$) W 2 $ + £>(/i£>$) 



-i(V$) 2 -^V-(VMV/>) 2 ), (2) 



where x= (x, y), V = (d x ,d y ), h is the surface displacement away from planarity, $ is the 
value of velocity potential at the free surface (the surface velocity potential), 7 = 2vk\jioo 
is a linear damping coefficient (7 ~ 5 2 <^ 1), G = g k /u>Q, T = Tk^/ (pui^), and / = 
<7zo&o/(4k>o) i s the driving amplitude, and V is a Fourier- integral operator, which is defined 
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for an arbitrary function u(x) by T>u(x.) = |k|u(k) exp(zk ■ x)<ik with u(k) the Fourier 
transform of u(x). Notice that Go + T = 1 by definition. Only those viscous terms that are 
linear in the surface variables are retained in Eqs. ([!]) and Nonlinear viscous terms, which 
were referred to as nonlinear dissipative effects in the introduction, have been neglected. 
Both the asymptotic analysis that follows and the numerical results presented later are 
based on Eqs. ([J) and (Q). 

Linear analysis indicates that for 7 1 the planar surface becomes unstable at f c = 7, 
and the critical wavenumber is k = 1. We seek nonlinear standing wave solutions near 
threshold (e = (/ — 7V7 <C 1) and expand Eqs. ([]]) and (0) in e 1 ^ 2 with multiple time 
scales, /i(x, t, T) = e l l 2 h x (x, t, T) + eh 2 + e^ 2 h 3 + • • • and $(x, t, T) = e 1 I 2 § 1 (x, t, T) + e<& 2 + 
e 3 / 2 $3 + • • •, where T = et. At O^e 1 ^ 2 ), we consider a set of N standing wave modes with 



critical wavevectors ±/cj; then the linear solution reads, 



hv 



f 



\ N 

cos t + ^ sin 3tj E 



Aj(T) exp(ifc 3 --x) + c.c. 



(3) 
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$1 = I — sin t + / cos t + cos 3t 



JV 

x E exp(i/cj-x) + c.c. 



(4) 



where the complex amplitudes Aj are assumed to vary in the slow time scale T. Notice 
that the linear solution contains not only the subharmonic responses but also the terms that 
are proportional to the driving force /. These latter terms arise because / is finite for the 
expansion in e. Terms that are of higher order in / have been neglected. In what follows, 
we shall replace / in the linear solution by 7 since / = 7(1 + e). 

At 0(e), there is no solvability condition; however there are resonant interactions that 
have to be taken into account. A particular solution for h 2 can be written as 



N r 

h 2 = J2\ H itt) [ A J A i ex P ( l C k j + h) ■ x) + 
3,1=1 L 



c.c. 



+H jt -i(t) [A,A\ exp (i(kj - fc,) 



c.c. 
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where the Hji(t) are unknown functions that satisfy, 



duHji + 2 7v /2(l + CjfidtHfl + [Go + 2r (l + c jt )] 

x y/2(l + Cj i)Hji = Fff cos 2t + FjP sin 2t + ■ ■ ■ , (6) 

where = cos% = fcj • ki, and and are proportional to the amplitudes AjAi. 
Equation (P) is the equation of an additively forced harmonic oscillator with friction. When 
the "natural" frequency of the "oscillator" equals to the driving frequency, resonance occurs. 



This condition reads, [G + 2T (1 + Cji)]y2(l + Cji) = 4. Due to the nonzero damping 

coefficient, this resonance results in a finite value for Hji that is inversely proportional to 

the damping coefficient. Since the RHS of Eq. (^) is proportional to AjAi, there are three 

waves involved in this resonance, namely, standing wave modes Aj and Ai, and mode B with 

wavevector kj + fcj. Therefore, Eq. (|6|) describes a three- wave resonant interaction. Note 

that the wavenumber for mode B is away from the critical wavenumber ko = 1; thus mode 

B is a linearly stable mode. The relevance of triad resonant interactions to pattern selection 

can be understood intuitively. Since the resonant growth of linearly damped mode B is 

at the expense of the reduction of the amplitudes Aj and Ai, the growth of modes Aj and 

Ai are less favored than any other mode A m that does not participate in a triad resonant 

interaction. In other words, the system tries to avoid critical modes that participate in 

triad resonant interactions. As shown below, triad resonant interactions strongly influence 

pattern selection through coefficients of cubic nonlinear terms in the amplitude equations 

in a way that is consistent with the above argument. Finally, for Tq < 1/3, triad resonance 

is not possible. For T = 1/3, wavevectors of the three resonating waves are in the same 

direction (8^ = 0). As r is further increased, off also increases. For purely capillary waves 

(r = 1), reaches the maximum value of 9$ ~ 74.9° or Cji = 2 1 / 3 — 1. 

At 0(e 3 / 2 ), we obtain the standing wave amplitude equations (SWAEs) from a non-trivial 

solvability condition, 

BA- r N 1 

y = -yA j - 7S(l)Kf + 7 E giCflW Aj, (7) 
where j = 1, 2, • • • , N, and g(l) and g(cji) are given in Ref. \T2\. There are two kinds of 



contributions to g(l) and g(cji). One is from the linear viscous terms in Eqs. ([I]) and (pj). 



The other is due to the parametric driving force and proportional to the driving amplitude 
/. These two kinds of contributions appear together in Eq. (|7|) since we have set / = 7 
at the linear order. The latter contribution is directly related to the terms proportional 
to / in the linear solution (Eqs. (0) and (^)), and it provides an amplitude-limiting effect. 
The nonlinear interactions between the primary subharmonic modes and terms related to 
the driving force in the linear solution produce terms that are out of phase by 7r/2 with 
the primary subharmonic mode, and thus can contribute to saturate the wave amplitude. 
An important point is that this amplitude-limiting effect results from the forcing term, but 
not from a dissipative term. As a result, this effect is also important even for Hamiltonian 
systems. 



It can be shown that g(l) > JT2|, which indicates the bifurcation to the standing 



wave state is supercritical. We rescale the amplitudes as yg(l)Aj and the coefficients as 
g(cji)/g(l). In what follows, we shall only refer to the scaled amplitudes and the scaled 
coefficients; but use the same notation for them as for the unsealed ones. Note that the 
scaled coefficient g(cji — > ±1) = 2. Figure 0(a) and (b) shows the scaled function g(cji) 
for two different values of the damping coefficient 7 and T = 1. The maxima in g(cji) 
around Cji = 0.26 (9ji = 74.9°) correspond to the triad resonance for purely capillary waves. 
The function g(cji) for capillary-gravity waves of Tq = 1/3 is shown in Fig. |Tj(c) and (d). 
Since the triad resonant interaction occurs among waves with their wavevectors in the same 
direction when r = 1/3, the resonant peaks (or maxima) are at Cji = ±1. Instead of large 
peaks at cji = ±1, we see that g(cji) has a wide flat region around Cji = and reaches very 
small positive values for small 7 due to scaling. In all cases, the effect of triad resonance is 
weaker for larger values of 7 as expected. 
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FIG. 1. The coefficient g(cji) as a function of Cji for purely capillary waves (Fo = 1) with the 
linear damping coefficient 7 = 0.02 (a), and 7 = 0.1 (b). The same coefficient for gravity-capillary 
waves of To = 1/3 with 7 = 0.02 (c) and 7 = 0.1 (d). 

The issue of pattern selection can be discussed by noting that Eq. (|7|) is of gradient form 
l/jd T Aj = —dT/dA*. Apart from the trivial solution of Aj = for j = 1, • • • , N, Eq. (^) 
has a family of stationary solutions differing in the total number of standing waves N for 
which Aj 7^ 0. By considering the case in which the magnitudes of all standing waves are 
the same, Eq. (|7]) has the following solutions, \Aj\ = \A\ = (l + Y^Liy^j) 9( c ji)) ■ The 
Lyapunov function for these solutions are, T = — ^l^l 2 — — + J2i=ia^j) 9{ c jl)\- For 
N = 1 (parallel roll solution), T\ = —\. For iV = 2, we have either square (C12 = 0) or 
rhombic [cyi 7^ 0) patterns with JF 2 = —1/(1 + gicyi))- If we consider only regular patterns, 
i.e. pattern structures for which the angle between any two adjacent wavevectors kj and 
k J+ i is the same and amounts to n/N, we have either hexagonal or triangular patterns 
for N = 3. Regular patterns for N > 4 are two-dimensional quasi-crystalline patterns 
(or quasi-patterns [[□[] )■ A quasi-pattern has long-range orientational order but no spatial 
periodicity, thus analogous to a quasi-crystal in solid state physics. Such patterns have been 
already observed in experiments of Faraday waves in systems driven by two carefully chosen 



frequency components fL3 |, but not in the single frequency case analyzed here. Fig. |2| shows 
the Lyapunov function F N as a function of 7 for iV = 1,2, 3, 4, 5, 6, 7, 8, and T = 1/3 and 1. 
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For T = 1, patterns of square symmetry (N = 2) have the lowest values of Tn for all values 
of 7 < 0.2, in agreement with experiments. However, patterns of different symmetries are 
favored in different ranges of 7 for To = 1/3. This can be understood qualitatively by noting 
that the self-interactions of the critical modes are less favored than interactions among them. 
Consequently, pattern structures with large N are favored. 
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FIG. 2. The values of the Lyapunov function Tm [N = 1, 2, 3, 4, 5, 6, 7, 8) as a function of 7 for 
purely capillary waves (ro = 1) in (a), and for capillary-gravity waves of To = 1/3 in (b). Part of 
(b) with small 7 is shown in (c). 

We next present results from numerical solutions of the QPEs (Eqs. (0) and to check 
the stability of the quasi-patterns described above. An analytic stability calculation is far too 
involved. We use a Fourier-Galerkin spectral method [|nj with periodic boundary conditions 
in a square domain. Time discretization is of second order. We use the trapezoidal scheme 
for the linear terms and the second order Adams-Bashforth scheme for nonlinear terms. The 
nonlinear terms are calculated by a pseudospectral method by using Fast Fourier Transforms. 
Most of the numerical studies are done for a system size of 647r x 64tt, or 32Ao X 32Ao with 
the critical wavelength X = 2ir/k = 2n, and a total number of 256 Fourier modes is used 
for each axis. We use a time step At = 0.04 — 0.1. Extensive numerical studies have been 
performed for different values of the three dimensionless parameters To, 7, and e in the 
QPEs. 

FIG. 3. Configuration of /i(x, t) at t = 32000 is shown in gray scale for 7 = 0.02, To = 1/3, and 
e = 0.1. 
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In the case of purely capillary waves (r = 1), Eqs. (P and (Q) are integrated for 
7 = 0.02, 0.05, 0.1, 0.2, and small values of e = 0.02, 0.05, 0.1, from a random initial condi- 
tion (gaussianly distributed with zero mean and a variance of 10 -4 ) for field h(x, t = 0), and 
zero values for $(x, t = 0). Asymptotically regular patterns of standing waves with square 
symmetry are obtained for all the above parameters; thus verifying the results of the asymp- 
totic analysis presented earlier, and again in agreement with experimental observation. For 
capillary-gravity waves at To = 1/3, the QPEs are integrated for three different values of 
7 = 0.02,0.05, and 0.1 at e = 0.05, with the same initial conditions. For 7 = 0.1, the long 
time configuration of /i(x, t) is of approximate square symmetry. For 7 = 0.05, the long time 
configuration of /i(x, t) is of approximate hexagonal symmetry. For 7 = 0.02, and to avoid 
finite size effects (expected to be stronger for smaller damping coefficient 7 [jl2| , p~3[1 ) , we have 
performed numerical studies for a system of size 64Ao x 64Ao with 512 Fourier modes for each 
axis and a time step At = 0.1. The initial condition for h is a set of gaussianly distributed 
random numbers with zero mean and a variance of 10 -6 , and $ is set to zero initially. The 
configuration of /i(x, t) at t = 32000 is shown in Fig. ^. The structure factor of this config- 
uration has been computed and has eight peaks, which correspond to eightfold symmetry 
of the standing wave pattern. This can also be seen by viewing Fig. |3| at a glancing angle. 
Thus, the basic prediction of the SWAEs is confirmed. Finally, we remark that the possible 
confirmation of stable pattern structures of even lower symmetry that are predicted by the 
SWAEs for smaller values of the damping coefficient 7 would require systems of larger size. 

This work was supported by the Microgravity Science and Application Division of the 
NASA under Contract No. NAG3-1284. This work was also supported in part by the 
Supercomputer Computations Research Institute, which is partially funded by the U.S. 
Department of Energy, contract No. DE-FC05-85ER25000. 
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